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ABSTRACT 

We study Minkowski Functionals as probes of primordial non-Gaussianity in the Cos- 
mic Microwave Background, specifically for the estimate of the primordial 'local' bi- 
spectrum parameter /^^ i with instrumental parameters which should be appropriate 
for the Planck experiment. We use a maximum likelihood approach, which we couple 
with various filtering methods and test thoroughly for convergence. We included the 
effect of inhomogeneous noise as well as astrophysical biases induced by point sources 
and by the contamination from the Galaxy. We find that, when Wiener filtered maps 
are used (rather than simply smoothed with Gaussian), the expected error on the 
measurement of f^^ should be as small as A/^^ — 10 when combining the 3 channels 
at 100, 143 and 217 GHz in the Planck extended mission setup. This result is fairly 
insensitive to the non homogeneous nature of the noise, at least for realistic hit-maps 
expected from Planck. We then estimate the bias induced on the measurement of 
/nl by point sources in those 3 channels. With the appropriate masking of the bright 
sources, this bias can be reduced to a negligible level in the 100 and 143 GHz channels. 
It remains significant in the 217 GHz channel, but can be corrected for. The galactic 
foreground biases are quite important and present a complex dependence on sky cov- 
erage: making them negligible will depend strongly on the quality of the component 
separation methods. 

Key words: Cosmology: Cosmic Microwave Background, Minkowski functionals; 
Methods: statistical, numerical 



1 INTRODUCTION 

The Planck satell it e was launched in Ma y 2009 
jTauber et al.l I2OI0I : iPlanck Collaboration et~aLl boilal ) 
to observe the microwave sky and particularly the Cosmic 
Microwave Background (hereafter CMB). It will provide 
major pieces of information about the early evolution of the 
universe and the origin of structures thanks to its unprece- 
dented combination of resolution, sensitivity, and spectral 
coverage. One of the important objectives of this mission is 
to bring new constraints on inflation theory and primordial 
non Gaussianities. In the simplest inflationary models based 
on a single slowly rolling scalar fleld, prim ordial fluctuations 
should be only we akly non Gaussian iMaldacenal l2003l : 
iBartolo et all I2OO6I '). Yet in more general models a much 
higher level of non-Gaussianity (hereafter NG) is expected. 
We can cite models with multiple scalar fields, features 
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in the infiaton potential, non-adiabatic fluctuations, non- 
standard kinetic terms, warm in flation, deviat i ons f rom 
Bunch-Davies vacuum (review in iBartolo et al. liooi) or 
completely different mechanis ms, for example topolog ical 
defects such as cosmic strings l|Kaiser fc Stebbins|[l984l ). 

Non-Gaussianity is often parametrised in a phe- 
nomenological way by the non-linear 'local' coupling pa- 



rameters 



which appear in the perturba- 



tive development of the primordial curvature p erturbation 
(|Komatsu fc Spergel|[200ll : IOkamoto fc Hull2002l '). 

^ix)^Mx) + f^d^Ux)-{c^Ux)))+g^^<t>Ux) + ..., (1) 

where 4il{x) is the linear Gaussian part of the Bardeen cur- 
vature. We will focus here on constraining the first parame- 

ter, 

The latest results from WMAP (|Komatsu et al.l 1201 ll ) 
with 7 years of data provide the constraint 

= 32 ± 21 (2) 

at the 68% confidence level. With Planck, we expect the 
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at best of the order of 2.5 from 
both temperature an d polarisation l|Yadav fc Wandef3l201(]| : 
ISefusatti et~all l2009h or of the order o f 5 for the temper- 
ature only (|Komatsu fc SpergelHiooH ). These constraints 
on were obtained with bi-spectrum m easurements 

ijKomatsu et al.ll2005l : IVadav et al.lfioOTl . l2008l l. The CMB 
bi-spectrum is the harmonic transform of the three-point 
correlation function and it has been shown that this is the- 
oretically an optimal estimator for /^^ as it satur ates the 
Cram er-Rao inequality for a weak non Gaussianity (jBabichl 

[iooi). 

However, alternative statistics to the bi-spectrum have 
been also developed, that can at least serve as checks and 
diagnoses of the results obtained from the bi-spectrum. 
One of the main reasons to study various probes of non 
Gaussianity is indeed that they are affected differently 
by different systematics and contaminants such as inho- 
mogeneous noise and foreground residuals induced by our 
galax y and point sour ces, as well as secondary anisotropies 
(e.g., lAghanim et al.l [2008, and references therein) such 
as integrated Sachs- Wolfe e ffect, Sunyaev-Zeldovich effect 
l|Sunvaev fc ZeldovichI IToT^ ) and weak gravitational lens- 
ing, which can all contribute to non Gaussianity of the 
CMB in a non trivial way. In particular, the most seri- 
ous bias on f^^ from secondary anisotropies arises from 
the coupling between weak gravitatio nal lensing and in- 
tegrated Sachs-Wolfe effect (see iGoldberg fc Spergell 19991: 
Serra fc Cooravll2008l: iHanson et al.ll2009l : iMangilh fc Verdd 



20091 : iMunshi fc Heavens! l2010l l. which adds complexity to 



the analyses. 

Here we focus on a complete set of topological tools, 
the Minkowsk i Functionals (here after MPs) introduced in 
cosmology by iMecke et al.l (|l994l ) . MPs describe the mor- 
phological features of random fields over excursion sets, 
i.e. regions where the field exceeds some threshold level i'. 
These well-known probes of primordial non-Gaussianities 
JSchmalzing fc Buch ertl Il997|; IWinitzki fc Kosowskvl 1 19981 : 
ISchmalzing fc Gors ki 1998) have been widely used in two 
and three dimensions, for instance on WMAP CMB data 



Hikage et al.ii200Ci, 2008|) and on the SDSS galaxy catalogue 



Park et al.ll2005l : iHikage et al.ll2006l l 



For CMB studies, MPs provide a nice complement to 
the bi-spectrum for several reasons. Firstly, at variance with 
the bi-spectrum, which evolves in harmonic (or Fourier) 
space, MPs are defined in real space, which makes a robust 
implementation for MPs in practice much easier than for the 
bi-spectrum. Secondly, MPs are sensitive to the full hierar- 
chy of higher-order correlations, instead of third order only, 
and can provide additional information on all the non-linear 
coupling parameters beyond the sole /nl. Because they are 
measured on excursions of the density field smoothed with 
an isotropic window, MPs only probe angular averages of 
higher order statistics, leaving out the angular dependences, 
at variance with the bi-spectrum. This explains why MPs 
are less optimal than the bi-spectrum in terms of disentan- 
gling between various models of weak NG. However, the very 
nature of these angular averages eases considerably statis- 
tical analysis and reduces the number of parameters while 
performing maximum likelihood analysis. As a result, even 
though suboptimal, MPs represent powerful tools of investi- 
gation of NG. For instance, the first limit on the primordial 
non-Gaussianity in the isocurvature perturbation was mea- 



sured with MPs, even though theoretical predictio ns where 
computed on the bi-spectrum (jHikage et al.l 120091 ) . Note fi- 
nally that the constraint on /nl obtained from the analysis 
of WMAP3 with MPs reads 



11 ±40 



(3) 



at the 68% confidence level (|Hikage et al.|[200^ . a rather 
competitive result, after all, compared to the bi-spectrum 
constraints from WMAP7 (eq. [2} . 

In this paper, we investigate constraints on f^^ ob- 
tained with MPs from a practical point of view, in order to 
ease at best future application of MPs to Planck data. After 
detailing the MPs estimators and the Bayesian method we 
implemented, we review the effects of the main systematics 
that can affect/bias the results. These systematics can be 
of instrumental nature -inhomogeneous noise, beams- or of 
astrophysical nature: the observed microwave sky contains 
foreground emission from the Milky Way and from extra- 
galactic sources. 

Our Galaxy is indeed a strong source of contamination. 
It is generally accounted for by masking and/or by using 
component separation methods and by assuming that the 
final residual bias on / is negligible comp ared to error 
bars (jHikage et al.ll200^ : lK"omatsu et aT1l2003h . Nevertheless 
it has been shown that component separati on leave Galactic 
features in CMB maps dChiang et aLlbOOSlV Indeed the var- 
ious methods of component separation ( Leach et al.ll2008h 
need to be controlled at the new level of accuracy on f^-^ 
expected for Planck. Another approach is to use fore ground 
templates and to marginalise over them (Koinats u et al.l 
|2002| ,[ 201 ll ). Here, we chose to analyse the behaviour of MPs 
when using masks and a naive model of component separa- 
tion. For that, we used foreground templates provided by the 
Planck Sky model (PSM) code (Delabrouille et al, 2012). 

Point sources are mostly unresolved galaxies (at least in 
low-density regions of Galactic emission), some emitting a 
signal sufficiently high to be detected individually, the others 
forming a diffuse unresolved background. Point sources con- 
sist first of radio-galaxies, active galactic nuclei which emit 
in radio frequencies with synchrotron process. They can also 
be dusty starburst galaxies which are observed via the ther- 
mal emission of their dust heated by the ultra-violet emis- 
sion of young stars. Current instruments are not able to de- 
tect all these galaxies individually and the integrated emis- 
sion of all the faint galaxies form a diffuse background, the 
Cosmic Infrared Backg round (GIB) which has been recently 
studi ed and observed (jPuget et ah 19961 : iLagache fc Pugej 
2m^, iPlanck Collaboration et al.ll2011? ). The strong point 
sources are accounted for by masking them but the fainter 
ones can still induce biases in NG studies. T he point source 
bi-sp ectrum has already been measured (iKomatsu et al.l 
[201^ and its effect in the est imates of evaluated, 
for bi-spectrum measurements llBabich fc Pierp aoli 20o2; 
ISerra fc Cooravll2008l : llacasa et al.ll2012l ). 

In addition to these contributions, as mentioned above, 
CMB contains secondary anisotropies imprinted between 
the surface of last scattering and present time. While we 
leave the treatment of these secondary anisotropies for fu- 
ture work, we shall treat in detail the contamination from 
our galaxy and point sources. 

This paper is organised as follows. In section 2, we re- 
view the Bayesian method we use to optimise the constraint 
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on /j^j^ from MFs. In section 3, we study the effects of pixeli- 
sation and filtering. Section 4 deals with the effect of inho- 
mogeneous noise. In sections 5 and 6, we study astrophys- 
ical sources of systematic effects: first point sources, then 
Galactic foregrounds. Section 7 summarises and discusses 
the results. For completeness, some additional technical de- 
tails can be found in Appendix A, which discusses analytic 
predictions for MFs in the weakly non Gaussian regime. Ap- 
pendix B, which details the algorithm used to measure MFs, 
Appendix [Cj with a convergence study of our x^, and Ap- 
pendix D, which provides some details on the Planck simu- 
lation experiments performed in this work. 

The pixelisation scheme adopted in this paper is the 
same as in Planc k processed maps, namely HEALPijQ 
ijGorski et al.ll2005h . 

Note finally that there are a number of tables to illus- 
trate the results. Some of these tables are purposely incom- 
plete to lighten the presentation. 



as Be tti numbers or as cold/hot spots. Ichingangbam et al] 
l2012h . Analytic formulae (theoretical expectation values) for 
the quantities Vk are known for Gaussian and weakly non 
Gaussian fields and are summarised in Appendix [X] An in- 
teresting property of these functionals is that the pure spec- 
tral dependence can be factorised out: 

V^:{u) = A^Vi:{u) (4) 

where the amplitude Ak is determined by the shape of the 
power-spectrum of the field fluctuations. The renormalised 
functionals Vk depend only on the non Gaussian corrections, 
i.e. on the behaviour of correlations functions beyond second 
order. The analysis of allows us to focus on non Gaus- 
sianity, which is the goal of this paper. 

In what follows, we denote by X an estimator of the 
quantity X. We calculate Vk = Vk/Ak on a pixelised map 
by proceeding as described in Appendix [Bl to estimate Vk- 
The quantity Ak itself is obtained by direct measurement of 
(To and (Ti on each map (see Appendix [K\ for definitions). 



2 METHOD 

In this section, we set notations by first recalling basics 
of Minkowski functionals and their measurement inside an 
excursion of varying height (§ I2.I|I . Then we detail the 
Bayesian approach we use to estimate f^^^ I2.2|l . by com- 
paring the measurement y of a functional or a combination of 
functionals to its expected value. Within a set of reasonable 
simplifying assumptions, in particular weak non Gaussian- 
ity, this approach reduces to a simple test I2.2.1|l . The 
numerical estimate of this requires the accurate calcula- 
tion of the covariance matrix of our estimator y from a large 
set of Gaussian simulations, as discussed in § 12.2.21 Addi- 
tional details concerning the Gaussian simulations set up as 
well as the method used to test the convergence of the are 
given in Appendix[C] Once we have a numerically robust x^ 
test, we study its sensitivity to important parameters of the 
problem, such as the number of bins used to explore the ex- 
cursion and the range of the excursion levels, as well as the 
choice of the functionals and/or their combination (ij 12. 3|) . 



2.1 Minkowski Functionals 

For a two-dimensional field 5 of zero mean and of variance 
aoi defined on the sphere, and smoothed with a window of 
typical size £ (to be defined later), we consider an excur- 
sion set of height v = S/ao, i.e. the set of points where the 
field exceeds the threshold i^. In what follows, we shall study 
the topological properties of the excursion with four quan- 
tities, denoted by 14 {k = 0,1,2,3). The first three ones 
correspond to MFs: Vo is the fractional Area of the regions 
above the threshold, Vi is the Perimeter of these regions and 
V2 is the Genus, defined as the total number of connected 
components of the excursion above the threshold minus the 
total number of connected components under the threshold. 
The fourth one, V3 that will be called the Number of clus- 
ters, also noted A^dustcrs, is just the number of connected 
regions above the threshold for positive thresholds and re- 
versely for negative thresholds (these regions are also known 
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2.2 Testing non Gaussianity with a x^ 

To provide constraints on non- Gau ssianity, we adopt a 
Bayesian approach similar to that of iHikage et all l|2008h . 
comparing measurement of normalised functionals Vk on the 
data map under consideration to the "theoretical predic- 
tions" obtained from the average of measurements of Vk on 
a large number of non-Gaussian simulations, the non Gaus- 
sian part being simply proportional to f^^-^. We will focus 
here on the non-linear coupling parameter local f^^^ but the 
method can be applied to other types of non Gaussianities. 

To perform the measurements, we consider an ensemble 
of nbins values of Ui ranging from — fmax and -fi^max, defin- 
ing a vector Vk = {vk{vi), - ■ ■ , i'fc(j'nbi„a)} for each of the 
functionals. The statistics under consideration will then be 
a vector of 71 nfunctionais X ribins elements, y = Cfc if only 
one functional is used in the analysis, or y = {vi,Vj ,■■■} ii 
a combination of nfunctionais > 1 functionals is considered. 

The Bayes formula writes 



(5) 



In what follows, we shall take a flat prior for 
/pjL, with Pifjsii^) a constant, while the evidence 
/ P{y l/NL)^(/NL)'i./NL ^ill .jrist considered as a normal- 
isation. 



2.2.1 Construction of a x^ tsst 

We assume that the likelihood P(y\f^i^) is a Gaussian, 
which allows us to define a simple x^ test for J^^- This 
approximation is expected to be good in the regime where 
the fluctuations of each component yi of the vector y remain 
small compared to its ensemble average {yi), when one con- 
siders a large number of realisations of the random field. 
Keeping these fiuctuations small sets practical constraints 
on (i) the value of I'max, which should not be too extreme 
to avoid sensitivity to rare events, (ii) the smoothing scale 
I, that should be small enough to probe a sufficient num- 
ber of independent harmonic modes on the sky. Because we 
consider only small values of I in most of what follows our 
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main concern is the choice of fmax- We s haU restrict to th e 
range i^max < 5 (see also the discussion in lGott et al.lll990l ). 
With all these simplifications, the posterior becomes 



'P(/nlI y) oc exp 



X^(y, /n 



with 



x^(y,/NL) 



(6) 



(7) 



= E E ly^ - M.)] [yj - yjif.J] (8) 

where J/(/nl) i*' the model under test. It can be de- 
rived analytically i n the weakly non Gaussian regime (see 
iHikage et al.l |2008| . and appendix \^ or by direct mea- 
surements on simulations. Since we want to take into ac- 
count complex additional contributions such as inhomoge- 
neous noise and masks, we use here the more flexible sec- 
ond approach with non- Gaussian simulations provided by 
lElsner fc Wandelj (|2009l '). So ^(/nl) = (2/(/nl)) is the mean 
of y measured over m^^ maps with a level non-Gaussianity 
. In practice, we shall take 



mNG = 200. 



(9) 



Since non-Gaussianity is weak ijKomatsu et al.ll201lh . 
one can compute the covariance matrix C in a Monte-Carlo 
fashion relying on simulations of Gaussian maps of the CMB 
(to which we add, when relevant, beam effects, noise, masks). 
In that case, the posterior probability distribution function 
P(/pjL I y) is expected to be very close to a Gaussian, as illus- 
trated by Fig. [T]for the Perimeter. Indeed, in the weakly non 
Gaussian regime, yi{fj^i^) is, at first order, a linear function 
of /nl (see, e.g.. Appendix 

To estimate the most likely posterior value of /nl and 
an uncertainty, instead of trying to find iteratively the max- 
imum of -P(/nl) followed by an estimate of the local cur- 
vature using e.g. the Fisher information matrix to estimate 
the error, we prefer, for simplicity and robustness, to com- 
pute directly the average and the variance of the posterior 
distribution, -P(/nl)- This latter is estimated using a num- 
ber of equally spaced bins, fi, i — 1, ■ ■ ■ ,n^^, with, for the 
purpose of this paper n^^ — 21. From this approximation 
of the posterior, the average and the variance are directly 
estimated numerically via the simple formulae 



In 



j:.p 



(10) 

(11) 



where Pi is proportional to -P(/nl ~ .A)- 

Finally, one might consider performing a number mtost 
of realisations of the data y -also called in this paper "test" 
maps- to have a more accurate description of the "typical" 
expected posterior probability. However such a forecasting is 
not obvious, as we have to define a frequentist average over 
Bayesian jjuantities, to predict the typical value expected for 

and A/j^^ . Our choice, equivalent to Fisher information 
forecast in the pure Gaussian case, consists in estimating 
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Figure 1. Example of measured posterior probability for /j^^ 
when using the Perimeter on one test map (red points). The green 
curve corresponds to a Gaussian of average fj^^ and of variance 
Af^^, with /jjj^ and Af^^ given by eqs. HlOj l and lllip . For this 
particular example, we used ribins = 26 and t'max = 3.5. The 
actual value of fj^^ is fj^^ = 0. 



/pjL and Afj^j^ from the following averages 

mtcst 



test 

I— 1 

1 ^2 

it L . ^ 



mtcst 



(12) 
(13) 



where f^^,i and A/^^ ■ are obtained from the analysis 
of "data" map number i using eqs. (|10p and (|ll|l . In what 
follows, we shall always take (for each individual value of 
/nl considered in the tests maps) 



mtcst = 200. 



(14) 



2.2.2 Convergence issues 

The covariance matrix is estimated from the average over m 
simulations, on each of which a y'^ vector of functionals is 
measured: 



Cij 



y, 



0( 



G 



-G 



(15) 



withyf = {yf). 

The calculation of the requires the inversion of the 
covariance matrix, an uneasy task, because C can be ill- 
conditioned, especially if m is not large enough. The num- 
ber of simulations required to have a good estimate of the 
X^{y,f^i^) function indeed depends on the number of bins, 
ribins, on the functionals under consideration, their number, 
rifunctionais, and ou the choice of i^max. The way we estimate 
m is exposed in detail in Appendix [C] The results of our 
analyses are summarised in table [T] where the minimum 
values of m required for having a better than two percent 
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^max — 5 




^max — 2 


'^bins 


= 11 


m 


Ai/ = 1 
= 7.0 X 10^ 


Au = 0.7 
m = 8.0 X 10^ 


Au = 0.4 
m = 8.0 X 10^ 


^bins 


= 26 


m 


Au = 0.4 
= 8.0 X 10^ 


Au = 0.28 
m = 9.0 X 10^ 


Au = 0.16 
m = 9.0 X 10^ 


^bins 


= 51 


m 


Av = 0.2 
= 9.0 X 10^ 


Au = 0.14 
m = 10.0 X 10^ 


Au = 0.08 
m = 10.0 X 10^ 


'^bins 


= 101 


m 


Au = 0.1 
= 12.0 X 103 


Au = 0.07 
m = 14.0 X 10^ 


Au = 0.04 
m = 15.0 X 10^ 



Table 1. Number of Gaussian maps m needed for a good con- 
vergence of the 'x^ statistic with an accuracy better than 2%. 
The calculations are performed here for the Perimeter, y = vi, 
and for various values of the number of bins n = n^ins a-nd of 
Umax, but these results would not change much for other func- 
tionals or combinations of functionals, as discussed in the main 
text. The details on the simulations are given in Appendix[C] For 
Umax = 5 we have removed extreme thresholds and reduced Wbins 
to {9, 22, 45, 89} respectively because the distribution of errors 
was not Gaussian for these bins which are too sensitive to rare 
events. 

accuracy on the obtained from the perimeter are dis- 
played for various realistic set-ups in terms of the number 
of bins and of the excursion range. Similar orders of magni- 
tude are expected for other functionals or when combining a 
set of functionals, as explicitly checked for the combination 
V1 + V2, ribins ~ 26 and t'max = 3.5, where the corresponding 
value found in Table [1] remains unchanged within 5%. This 
latter property comes from the fact that in practice, two 
different functionals are much less correlated than the same 
functional estimated at two successive bins of the excursion: 
combining functionals doubles the dimensions of matrix C 
but does not make it significantly more degenerate. In the 
subsequent calculations performed in this paper, we shall 
take 

m = 10 000, (16) 

a number of simulations sufficient for a good estimate of the 
for Jibin < 50 and Umax > 3.5. 

2.3 Sensitivity of the estimator 

In this section we discuss the sensitivity of our estimator, 
that is the uncertainty on the estimated /^^ 1 with respect to 
the number of bins, nbins, the excursion range, I'max, and the 
set of functionals used, whether it is the Area, the Perime- 
ter, the Genus, the Number of clusters, or a arbitrary com- 
bination of any of them. Our simulation setup is the same 
as in the previous paragraph (and detailed in Appendix [C)) . 
We checked, in this full sky configuration with homogeneous 
noise, that our estimator is unbiased, (/^l) — /nl (within 
numerical limits set by the finiteness of tting and mtcst). 

The results of our analyses are displayed in Tables [2] 
andO They show that the combination 

ax 7 ^bins )^ (3.5, 26) (17) 

is fairly optimal and shall represent our choice in the sub- 
sequent analyses. Note that it is important to notice that 





,1/2 


n = 11 


n = 26 


n = 51 


n = 101 


t^max — 


5 


29 


25.5 


25 


25 


i^max — 


3.5 


27 


25 


25 


25 


t'max — 


2 


37.5 


37 


37 


36.5 



Table 2. Investigation of the best combination for (i^max, ^bins)- 
This table gives the uncertainly in /j^^ for different values of the 
number of bins n = nbins a-nd the threshold fmax. The calcula- 
tions are performed for the Perimeter, which is the most sensitive 
to /j^j^ . The results would be nearly the same for the Genus and 
^clusters I while the Area is quite insensitive to n and fmax for 
the range of values tested here. From this table it is fairly easy to 
conclude that the combination (i^max, "bins) = (3.5, 26) is close to 
optimal. In particular it is not necessary in practice to go beyond 
f^bins ^ 26. Note for completeness that the error bars computed 
here correspond to the specific set up of Appendix \C\ (see also 
caption of next table). 



Functionals Vj. 


^ NL ' 


Vb 




44 


Vi 




25 


V2 




25 


V3 




28 


V1 + V2 




20 


V1 + V2- 




19 


V1 + V2- 


h V3 -f Vb 


18.5 



Table 3. Sensitivity of the estimator: error bars on the measure- 
ment of /jjL lor each functional, and for various combinations. 
The calculations are performed assuming n = ifynctionals x "bins 
with ribins = 26 and i^max = 3.5. Note for completeness that 
the error bars computed here correspond to the specific set up of 
Appendix [C] (see also Appendix [DJ : extended mission for the 3 
channels and noise filtering with a Gaussian window of full width 
at half maximum Sp^jjj^j = 10' and HEALPix resolution param- 
eter A^sidc = 2048. 



assuming t he covariance matrix t o be diagonal, as d one for 
instance in lEriksen et al] ^2004 ): iGott et al.1 (|l990l ) ts not 
a good approximation in the case of /„l and decreases by 
more than a factor two the constraining power of the x^- 

The comparison between various functionals lead to the 
following ranking, in term of discriminative power: 

Perimeter > Genus > Afdustcrs 3> Area. (18) 

While the Perimeter, the Genus and the Number of cluster 
present the same order of sensitivity, the Area is about twice 
less discriminative than them. Most of the information on 
/„L can be extracted by a combined analysis of the perime- 
ter and the genus, Vi + V2, with a little improvement when 
taking into account the number of clusters, Vi + V2 + V3, 
but the area does not carry significant pieces of additional 
information compared to the three others. 
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3 FILTERING AND SMOOTHING 

The measurement of Minkowksi Functionals requires 
smoothing of temperature maps in order to remove the con- 
tribution of the noise. This smoothing is performed at var- 
ious scales on a pixelised map in order to extract all the 
statistical information available. First, one has to deal with 
discreteness effects brought by the pixelisation. In § 13.11 we 
show that for practical measurements, it is not needed to 
have a large value of IfEALPix resolution parameter A^'sidc 
to extract all the relevant information even if this means a 
pixel size comparable to the size of the smoothing window: 
this is because we explicitly account for these pixelisations 
effects in the model. In § 13.21 we consider Gaussian smooth- 
ing and show that most information on /^^ obtained from 
MFs is present at the smallest possible scales, i.e. at scales 
comparable to the size of the beam of the instrument. Note 
however that this result stands for local f^-^ and might vary 
for different types of non Gaussianity. Finally, § 13.31 deals 
with Wiener filtering for the field, its first and second deriva- 
tives. We show that the results obtained with a simple com- 
bination of Wiener filters set much better constraints on /j^^ 
than naive Gaussian smoothing at various scales. 

In what follows, smoothing scale will be expressed in 
units of the full width at half maximum (FWHM), Ol^^^ 
(see Appendix ID2I for more details). 

3.1 Pixelisation effects: choice of A'^sidc 

In principle, the pixel size should be small compared to the 
smoothing kernel size, in order to avoid systematic errors 
intr oduced by the discrete nature of the pix elisation (see, 
e.g. IColombi et"al1l2000l : iNovikov et all 120061 '). We checked 
in practice that no bias is introduced on the measurement 
of f-^-^^ even when the pixel size becomes comparable to the 
smoothing window size, because the defects of the pixelisa- 
tion are present as well in the model. However, introducing 
larger pixels is similar to coarse graining and decreases the 
effective level of statistical richness, which in turn increases 
the uncertainty on f^-^ . Furthermore, a too large pixel size 
would simply make the additional filtering operation inop- 
erative: instead, what we would get in that case is simply 
the dominant part of filtering to be a convolution with a 
top-hat of the pixel shape. 

Table 2] gives the error obtained on f^-^ when using 
the combination of all four functionals, for various values 
of HEALPix resolution parameter A'sido, different levels of 
noise and a Gaussian smoothing with = 5'. For 

r W rl iVl 

Planck purpose, we also find that the improvement brought 
by the A'side = 2048 resolution is tiny, as expected for higher 
level of noise, and this table shows that it is not needed in 
practice to go beyond A'sidc = 1024, which is rather handy 
computationally speaking. From now on, unless specified 
otherwise, we shall assume 

iV.idc = 1024 (19) 

for all subsequent analyses. 

3.2 Gaussian Smoothing 

Smoothing with a Gaussian kernel depends only on angular 
scale, 6'p^jj„. A priori, each type non-Gaussianity is char- 



Noise 


A^sidc = 512 


Af.idc = 1024 


A^sidc = 


2048 




^max = 1000 


fmax = 2000 


^max — 


3500 


0.33 /xK.deg 


19.5 


16.5 


16 




0.5 /xK.deg 


20 


16.5 


16.5 




0.7 ^K.deg 


20 


17 


16.5 




1.1 /^K.deg 


20 


18.5 


18 





Table 4. Sensitivity (A/^^)-*^/^ versus A^sido for various levels of 
noise (see Appendix ID)| . when using the combination Vb-I-Vi-|- V2 + 
V3 of all four functionals and Gaussian smoothing with = 
5'. Other parameters used to perform the simulations are n^ins = 
26, Vmax = 3.5, m = 10 000 Gaussian maps for the covariance 
matrix, m^Q = 200 reference maps with different levels of f-^-^^ 
for the model and 200 test maps with /j^^ = 0- Note that, in 
the weakly non Gaussian regime considered here, the forecast 
(A/j^^)^/^ does not depend significantly on the actual value of 
/fjL ■ For each value of Afgide , calculations in harmonic space are 
band limited to i ^ ^max — 27Vsidc. This cut-off at ^max does not 
affect significantly the results. 



acterised by a specific scale range where the sensitivity of 
the estimator of f^-^ is the best and we have to stay aware 
that we are limited here to only one particular case of non 
Gaussianity, although quite typical. 

We tested different smoothing scales, including even 
very small scales (^p„jj^j = 5'), smaller than size of the 
beam (fi'p^jj^ ~ 7.2' in the case of the combined channels), 
and just above the size of the pixel, which might look awk- 
ward, but in fact does not introduce any significant bias on 
the measurement of f^-^ , as already argued in § 13.11 when 
discussing about pixelisation effects. 

One issue about smoothing at large scales is that it 
reduces the number of independent modes available on the 
sky. This in turns reduces the quality of the measurement 
of the MFs in the tails (large values of and can make 
the likelihood function non Gaussian, as discussed in the 
beginning of § 12.2.11 particularly if i^max is too large. For 
instance, with i^max ~ 4, the Gaussian assumption for the 
likelihood is legitimate for O^^.^^ = 5' but not for O^^.^^^ = 
40'. Here we checked that our analysis was still valid and 
well converged at smoothing scales as large as ~ 
for our default choice for the parameters, i/max ~ 3.5, ribins ~ 
26 and m = 10 000. 

Table [5] summarises the results of our analyses for the 
case when the combination of all functionals is used at vari- 
ous scales or various combinations of scales. Note that when 
combining 4 scales, we have a large number of entries in 
the data vector, n = ribins x nfunctionais x riscaies — 400, but 
we checked that this did not affect the convergence of the 
calculation of the covariance matrix with m = 10 000. The 
results of Table[5]show that there is no significant statistical 
information for > 20'. Most of the signal is captured 

by the combination fp^^^ = 5' with 0^^^^ = 10'. 

Note that the scales we consi der here are much sm aller 
than those considered in Table 3 of lHikage^ et al.1 (|2006! i that 
correspond to dp„^^ = {11.75', 23.5', 47'}, which explains 
the better constraints we obtain on f^^ . 
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FWHM 



5' 
10' 
20' 
40' 

5' 4 
5' 4 



10' 
10' 



20' 



5' + 10' + 20' +40' 



16.5 
20 

26.5 
34 
14 

13.5 

13.5 



Table 5. Sensitivity of the estimator of /nl versus Gaussian 
smoothinEf scale 9^ when the combination of all functionals, 

Vb + Vi + V2 + V3, is used. The forecasted quantity (A/j^^)^''^ 
is given for the combined 100, 143 and 217 GHz channels in the 
extended Planck mission (see appendix iDt . 



3.3 Wiener filters 

Now we consider Wiener filtering, which is an optimal way 
of recovering the signal of a data map D = S + N in case 
both the underlying signal S and the noise TV are Gaussian. 
The Wiener filter writes, in harmonic space. 



Wm = 



(20) 



where Ct is the forecasted power-spectrum of the signal we 
want to analyse and Ni is the known power spectrum of the 
noise. Here, signal will refer to the dominant, Gaussian, cos- 
mological part of the CMB. Obviously, this is sub-optimal, 
since we are after the non Gaussian part of the CMB, while 
all the rest should be considered as noise. However, using a 
Wiener filter designed for extracting the non Gaussian signal 
would require a stronger prior, making the measurement of 
/nl potentially more accurate but only for a restricted class 
of non Gaussianities. 

Similarly, one can define in harmonic space a "deriva- 
tive" Wiener operator. 



Ci 



(21) 



and a "second derivative" Wiener filter optimal for recover- 
ing the Laplacian of the signal, AS, 



Wd2 = £{£ + 1) 



(22) 



In practice, we use a smoothed version of the 
quantity Ce/{Ci + Ne), a "Wiener- like" function used 
for component separation method (CMB removal) in 
iPlanck HFI Core Team et al.l (|201lf ). Before applying the 
Wiener filters we also correct the map for the Gaussian beam 
corresponding to the channel configuration. The three re- 
sulting filters, Wm, Wdi and Wd2 are represented on Fig. [21 

Table [6] summarises the results of our analyses using 
them. It shows, not surprisingly, that Wm alone does much 
better than Gaussian smoothing (Table [SJ and a signifi- 
cant improvement is obtained when combining Wm with 
Wdi resulting in an overall reduction of the error on /nl 
of 30% compared to the best results obtained with Gaus- 
sian smoothing. On the other hand, Wd2 does not bring 
anything interesting, but this is somewhat expectable: for 




1000 1500 2000 2500 3000 

i 



Figure 2. The three Wiener filters, Wm, Wui and Wd2 given 
by eqs. ([21} and (|22j. 



a stationary and isotropic random field, there is no correla- 
tion between the field and its first derivatives, but there is 
a strong correlation between the field and its second deriva- 
tives. 



4 INHOMOGENEOUS NOISE 

Due to the scanning strategy and the orientations o f the 
different horns and bolometers l|Dupac fc TaubedbOOSl ). the 
distribution of the noise in the raw sky maps is correlated, 
anisotropic and inhomogeneous. In what follows, we treat 
effects of inhomogeneity and anisotropy of the noise, by 
relying on hitmaps ( generated with the software MADAM of 
iKeihanen et al.l 120051 ). in which each pixel value represents 
the number of times the pixel has been observed by the satel- 
lite. Our modelling of the noise in each pixel of the map i 
then reads 



noise (i) = aisotropic noise X A/'(0, 1) x 



(bitmap) 
bitmap (i) 



(23) 



Note thus that, since correlations are neglected, anisotropy 
of the noise is modelled only partly, through the anisotropy 
of the hit map. Even more realistic analyses would take into 
account of the effects of correlations in the noise, which will 
be addressed in future work. Figures [3] and [4] show a hit map 
and the corresponding noise map for the 143 GHz channel in 
a simulation of Planck signal ( which used the characteristics 
of the instrument described in lPlanck HFI Core Team et al.l 
|201li). 

To analyse the impact of various levels of realism in 
modelling of the noise on the determination of f^^^ , we con- 
sider two configurations: one where inhomogeneous noise is 
included in all the steps of the analysis, and one where 
only each of the mtcst = 200 test maps has its own realisa- 
tion of inhomogeneous noise - relying on the same hit majQ 



The one obtained in the 143 GHz channel, to be specific. 
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Functional Wm Wbi Wr>2 Wm + Wbi Wm + Wu2 Wm + Wui + Wd2 



Vo 


51 


Vi 


14 


V2 


21 


Vs 


20 


V1 + V2 


13 



V0 + V1 + V2 + V3 12.5 32 27 9 12.5 9 



Table 6. Sensitivity (Af^^ ) ^/-^ when Wiener filters are applied to the data maps prior to MFs measurements. 
It is calculated in the framework of the extended mission for the 3 channels. 



Hitmap@143GHz 




^^^^B ^^^H 10000 



Figure 3. A typical hit map in Galactic coordinates, obtained for 
iVgide = 2048. Pixel values quantify the number of observations 
of the pixel. Areas near the ecliptic poles are observed several 
times more frequently than regions of the sky near the ecliptic 
plane. The lowest values are about 25 while the highest are about 
20 000. 



Noise© 143GHz 




-6.8 ^^^^^^^^^ '^^^^^ G.e 



Figure 4. A typical map of inhomogeneous noise, using equation 
II23I I on the bitmap of Fig. |3] The map is shown in units of its 
standard deviation. 



but with different random seeds. We test here the null hy- 
pothesis = but we checked that the conclusions do not 
change significantly for other values of /j^^^ . All the simula- 
tions are performed as detailed in Appendix |D] for full sky 
surveys with no foregrounds, A^'sidc = 2048, and a Gaussian 
smoothing with ^p„jj„ = 5' to filter out the noise, which 
corresponds to a rather (almost the most) pessimistic case 
in terms of inhomogeneous noise. 



The results of our analyses are summarised in table [7] 
and table[8]for four levels of noise which are likely to bracket 
the actual sensitivity of Planck. To understand the results 
displayed in the tables. Figures [S] and [S] compare the ef- 
fect of neglecting inhomogeneous noise to the presence of a 
"true" on the functionals for the nominal mission in the 
143 GHz channel. 

The Area functional, Vo, seems fairly insensitive to the 
effect of the inhomogeneity of the noise, which in turns 
makes the determination of /^^ from the Area quite ro- 
bust to that respect. This is not very surprising: the area is 
proportional to the cumulated one point distribution func- 
tion (pdf ) . The presence of inhomogeneous noise locally in- 
duces a convolution of this pdf with a Gaussian of vary- 
ing width depending on the value of the number of hits in 
the map. The effect of this convolution is negligible when 
the r.m.s. anoisc of the local noise is small compared to the 
r.m.s. ao of the field under consideration Here this is the 
case: the additional Gaussian smoothing with S^,„„.. = 5' 

'~' r W rl JVL 

reduces the typical local rms of the smoothed noise map to 
fnoiso ^ 4 X 10~^ whatever the channel considered, to be 
compared to ao ~ 4 x 10~^. 

The examination of Fig. [S] shows that other function- 
als are rather sensitive to the effect of inhomogeneous noise, 
which is also expected. Indeed, we can guess that the pres- 
ence of inhomogeneous noise augments the contrasts be- 
tween cold and hot spots compared to the homogeneous 
case, hence building up a signal in Vi, V2 and V3. However, 
this signal is also contained in the parameter Ak in eq. Q 
which tends to compensate for the effect on Vk- Hence, the 
normalised functionals, Vk, appear to be less affected by in- 
homogeneous noise than the "raw" functionals, Vk, as illus- 
trated by Fig.[6l There is still some rather significant residual 
signal, at least for the small smoothing scale considered here. 
However, the parity of the black curves in the different pan- 
els of Fig. |6]is opposite to that of the curves corresponding 
to a true primordial /^^ (green and red curves): we do not 
expect in that case the presence of inhomogeneous noise to 
introduce any bias on the measurement of /nl . These simple 
statements are confirmed by the examination of Table [T] On 
the other hand, the presence of inhomogeneous noise makes 
the uncertainty on f^^^ slightly larger, increasingly with the 
average level of noise, as shown in Table [S] Fortunately, for 
the three combined channels in extended Planck mission. 



Note that this argument would be valid as well for a non Gaus- 
sian noise. 



Non Gaussianity and Minkowski Functionals 9 



the efTects of the inhomogeneity of the noise become nearly 
negUgible and can be fairly ignored, as we shall do from now 
on. 



5 FOREGROUNDS I: POINT SOURCES 

"Point sources" refer to the (large) number of radio and 
infra-red galaxies that are detectable at the CMB frequen- 
cies. These galaxies are in general not all resolved by CMB 
all-sky experiments like WMAP or Planck, even if the 
brightest objects are detected individually. The faint ones 
contribute to an inhomogeneous sky background. In this sec- 
tion, we test for the first time the effect of these sources on 



ous studies relying on MPs analyses (e.g 



mctionals. In previ- 
■ iHikage et allboosl : 



the estimate of /^^ with Minkowski Fu nctionals. In previ- 

iNatoli et al.l |2010|) , the effect of point sources was indeed 
supposed to be completely subtracted off by an appropriate 
masking or simply negligible compared to error bars. 

This section is divided into two parts: §[S7T] describes in 
details the simulation pipeline we used to perform our anal- 
yses while § 15.21 discusses the results of our investigations. 



5.1 Method 

To perform our simulations of the data, we need to com- 
pute accurately the contribution of point sources. To do 
so, we use the Planck S ky Model (PSM) code described in 
iDelabrouille et al.l (|2012t n. As reviewed in § IS.l.fl the point 
source contribution depends significantly on frequency. As a 
result, our simulations and analyses will consider separately 
the three cosmological channels, fOO, 143 and 217 GHz, in 
the extended mission configuration for the noise level (see 
Appendix [DJ. In particular each channel will have differ- 
ent masking treatment for the brightest point-sources. Since 
masks represent a crucial part of the treatment of point 
sources, we discuss about them in § 15.1.21 Other technical 
details about our simulations are provided in § 15.1.31 

5.1.1 Point sources simulations 



The Planck Sky Model (PSM) code ijDelabrouille et a l."2012') 
is specifically designed to simulate all relevant sky emissions 
at Planck frequencies, including secondaries and foreground 
emission, as they were known before the launch of Planck. 
In this paper we used only the part of the PSM that deals 
with point sources, the rest of the simulation pipeline being 
detailed in § 15.1.31 

Firstly, we use the PSM code to add radio sources, 
namely Active Galactic Nuclei (AGN), to our simulations. 
These AGNs are observed via their synchrotron emission. 
The PSM relies on numerous surveys of radio sources at 
frequencies ranging from 0.85 GHz to 4.85 GHz to model 
this emission. In regions not observed by surveys or with 
shallower observations, sources are copied from other re- 
gions, until a coverage down to at least 20 mjy at 5 GHz 
is achieved over the full sky. Then flux densities are ex- 
trapolated at all frequencies by using a power law ap- 
proximation for the spectra, of the form Su oc u"". For 

* |http : //www ■ ape ■ univ-parls7 ■ f r/ [ ^delabrou/PSM/ psm . html 



the spectral index a estimates, sources are classified into 
steep or flat spectrum class and a is drawn from a Gaus- 
sian distrib ution with mean and variance corresponding 
to i ts class (iRicci et al.l i2006) and matching WMAP data 
(Ben nett et all 120031 ) . in several frequency ranges. Besides, 



WMAP sources are accounted for in the simulations. Source 
counts at 5 and 20 GHz are found to be consistent with 
the model of 'Toffolatti et al. (Il99 8'l and an updated ver- 
sion of the model of d c Zotti et al.l (2005 ). These sources are 
known to contribute essentially at low frequencies, from 30 
to 90 GHz but t hey have been detected at higher fr equencies 
up to 217 GHz l|Planck Collaboration et al.li2011al lbh. 

We note that radio sources are nearly Poisson dis- 
tributed on the sky, so they essentially contribute to the 
CMB as a shot noise. 

Secondly, we add infra-red (IR) sources to the simula- 
tions. Indeed, at high frequencies (above 150 GHz), a ther- 
mal emission arising from dust heated by the UV emission of 
young stars also contributes. In addition to normal stars sur- 
rounded by a disk, numerous starburst galaxies which form 
stars at extreme rates contribute to this thermal emission. In 
the PSM, sour ces are taken from the I RAS Point Source Cat- 
alogue (PSC) (iBeichman et al. 1 19881 ) and the Faint Source 



Catalogue (FSC) (|Moshir et al. 



19921 ). The flux densities are 



extrapolated to Planck frequencies by adopting a model with 
modifled black body spectra and the gaps in sky coverage are 
flUed up using the same procedure as for the radio sources. 

Finally, we add to the simulations the Cosmic Infra-red 
Background (GIB) which is possibly the dominating compo- 
nent. Distant starburst galaxies are not all detected individ- 
ually and the cumulated emission from the fainter ones form 
a diffuse backg round of an i sotrop ies. The PSM adopts the 
count model of iLapi et al.l (|2006l ). which is consistent with 
SCUBA and M AMBO surveys . Sour ces are clustered follow- 
ing model 2 in lNegrello et al " (l2004h and the spat ial distri- 
bution follows the method of Gonzalez-Nuevo et al.l (|2005l ). 
then flux densities are extrapolated to all frequencies. The 
GIB power spectra of the simul ation agree sufficiently well 
with t hose measu red by Planc k (Planc k Collaboration et al.l 
l201lj ) and ACT (jPunklev et al. 2011) for our forecast stud- 
ies. 

The main difference between the distribution of radio 
and IR sources on the sky is that IR sources, being either 
observed as individual entities or as a background, are clus- 
tered in their host dark-matter halos. So their power spec- 
trum is not flat as for the radio sources. 



5.1.2 Masks 

Brightest point sources can be detected individually and 
can be masked properly when their flux density is be- 
yond a chosen detection threshold when compared to the 
level a — (Tnoisc of the underlying noise in the CMB data. 
Here, we create three sets of point source masks corre- 
sponding to three flux density cuts (referred simply as 
flux cuts). The choice of these different flux cuts has been 
mainly derived from the Early Release Compact Source 
Catalogue (ERGSC) published by t he Planck Collaboration 
(|Planck Collaboration et al.ll2011bl ). The flrst set of masks 
refers to sources beyond a lOo" level in the 3 bands 100, 
143 and 217 GHz, corresponding to respective flux density 
thresholds at 0.5, 0.33, 0.28 Jy as chosen in the ERGSC. 
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Noise (/^K.deg) 


1.1 




0.7 




0.5 


0.33 


Functional 


(/nl> (A/^)''" 


{/nl ) 




(/nl) 






Vb 




1 


38.6 








Vi 




-4 


23.6 








V2 




-0.4 


22 








V3 




-0.2 


24.7 








Vb + Vi + V2 + Va 


-4.4 18.8 


-3.8 


17.9 


-2.3 


17 


-0.8 16.7 



Table 7. Effect of neglecting the presence of inhomogeneous noise when estimating /j^^ s-iid the error bar on it in the 
null hypothesis, /j^j^ = 0- Each of the first two columns corresponds to a given level of noise expected in some specific 
Planck channel, respectively the 217 GHz and the 143 GHz channels, while the third and the fourth ones correspond to 
the combination of the 100, 143 and 217 GHz channels in the nominal and the extended mission case, respectively (see 
Appendix |D]| . The numbers in this table assume a Gaussian smoothing of the data maps with = 5'. 



Noise (/iK.deg) 




1.1 




0.7 




0.5 




0.33 


Functionals; V0 + V1+V2 + Vs 


(/nl> 


(A/^ )l/2 


(/nl) 


(A/^ )l/2 
^ NL ' 


{/nl ) 


(A/^ )l/2 
^ NL ' 


(/nl) 


(A/^ )l/2 
^ NL ' 


Configuration 1'^ 


0.19 


17.8 


-0.05 


16.5 


0.044 


16.2 


0.07 


16 


Configuration 2_^ 


-4.4 


18.8 


-3.8 


17.9 


-2.3 


17 


-0.8 


16.7 


Configuration 3_^ 


0.002 


25.5 


-0.44 


23 


-0.15 


19.9 


-0.15 


17 



Table 8. Effect of inhomogeneous noise when estimating /j^j^ s-^d the error bar A/j.^^^. Three settings are considered, as detailed 
below. The table is otherwise similar to Table [7] 



For the reference maps and the test maps, we introduce isotropic noise. 
^ For the reference maps, we introduce isotropic noise and for test maps we introduce inhomogeneous noise. 
For the reference maps and the test maps, we introduce inhomogeneous noise. 



The second set of masks concerns sources beyond the Scr 
level, which corresponds to the threshold choice in the clean- 
est parts of the sky of the ERCSC. The third set of masks 
corresponds to the 3cr level, which is not mentioned in the 
ERCSC, but we use it because we believe it represents a 
more appropriate set of masks for cosmological purposes. 
Indeed, in the ERCSC, the goal was to set flux density 
thresholds to have a sufRciently good signal to noise ratio 
for reliable analysis of the point sources properties. Our goal 
here is to remove the contribution from the point sources, 
which requires a much less stringent criterion on the quality 
of their detection. Furthermore, the ERCSC signal to noise 
level does not match that of the nominal mission and by no 
mean that of the extended mission. 



Each mask associated to an individual point source is 
a disk of radius 3 times the FWHM of the beam of the 
instrument in the channel under consideration. When adding 
up the contributions of all the sources, a certain fraction 
1 — /sky of the sky is masked, as illustrated by Fig. [T] The 
value of the sky fraction which is then used, /sky, ranges 
from /sky = 0.90 for the 100 GHz channel up to /sky = 0.99 
for the 217 GHz channel. These differences come from two 
factors on which the construction of masks depend: the beam 
width and the number of point sources detected beyond the 
threshold of interest. These two parameters decrease when 
passing from 100 to 217 GHz. 



PS@il43GHz 




0.0 ^^^^^K- ^^^H 1.0 



Figure 7. Masks of point sources at 3(7 (ERCSC refer- 
ence), designed for the 143 GHz map, drawn from the PSM 
llDelabrouille et al.ll2oT3) . 



5.1.3 Simulation pipeline 

To test the effects of point sources in the estimation of f-^-^ or 
more specifically the approximation of neglecting their pres- 
ence, we add their contribution only to the "data" (test) 
maps, y in eq. ((8]). We simulate mtest = 200 of these tests 
maps with /p™ = {-10,0,10,50}, where /p™ stands for 
the "primordial" /„l (to be contrasted later with other con- 
tributions to the effective f^^^ arising from biases induced 
by unaccounted point sources). The m = 10 000 Gaussian 
maps used to compute the covariance matrix C as well as 
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Figure 5. Effect of inhomogeneous noise on tlie "raw" functionals Vk for the nominal mission in the 143 GHz 
channel and a Gaussian smoothing with 6^„„„. = 5'. The relative difference between the Vi- measured in different 

^ FWHM ^ 

types of maps and the Gaussian limit is displayed as a function of v. Each panel corresponds to an individual 



functional. The curves represented on each panel are calculated by the average over 200 realisations, with /^^ 







and inhomogeneous noise for the black thick curve and with primordial non Gaussianity {f^^ 0) and homogeneous 
noise for the green and red curves, while Vp was computed alike with homogeneous noise. 



the rriNG = 200 non Gaussian maps used to calculate the 
model prediction j7(/nl ) in eq. (HI neglect this contribution, 
but have exactly the same treatment otherwise, including 
sky coverage and instrumental noise as detailed below. This 
way, our analysis will be able to confirm if appropriate mask- 
ing is enough to render the efTects of point sources negligible 
on the measurement of f^^ . 



by diffusive in-paintinglf| Then, the maps are smoothed with 
a Gaussian window of size Sp.j^jj„ = 5'. Finally, a galactic 
mask is applied to the maps, here corresponding to a valid 
fraction of the sky /sky ~ 0.80. The procedure to construct 
the galactic mask will be described in § [S] The complete 



The details of our simulation pipeline now follow. CMB 
maps are created first with the beam corresponding to each 
channel frequency / (see Appendix |D|, are supplemented 
with point sources (only for the test maps) convolved with 
the same beam and with the noise corresponding to each 
channel / for the extended mission. Next, point source 
masks are applied to the maps. These masks depend on the 
channel / - so the beam width &fwhm(/) is a parameter- 
and on the chosen flux cut /3. The punched holes are filled 



^ Choosing a lexical order (defined by HEALPix), the values in- 
side masked pixels are computed using the average over the values 
inside neighbouring pixels, when available, whether it is from an 
unmasked pixel or a pixel inside the mask that was calculated 
with the algorithm in a previous step. To achieve convergence, 
the process is reiterated a number nn of times. We take nu = 30, 
which is sufficient in practice for the mask size we have in our 
simulations. 
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Figure 6. Same as in Fig. (5] but for the normalised functionals wj., showing that these latter are much less sensitive 
to the effects of inhomogeneous noise than the "raw" functionals, Vjj. 



pipeline is summarised as follows: 

map = CMB *beam[6'pwHM(/)] 

+ foreground of sources * beam[0p\YHM(/)] 

+ noise (/) + point sources mask (/3, /) inpainted 

+ smoothing (Sp^^^ = 5' ) 

+ galactic mask (/aky). 



(24) 

We checked that the results derived in this section are qual- 
itatively the same for other values of d'^^^^^ and for the 
Wiener filters studied in ij 13.31 Of course, a quantitative cal- 
culation of the biases induced by point sources will require 
a new analyse each time a new filter is considered. 



5.2 Results 

Table |9] shows the estimate on f^-^ obtained for different 
channels as a function of source masking level and primor- 
dial non Gaussianity, /^l™- Again, a frequentist average of 
the posterior averages is performed over 200 test maps real- 
isations, and is noted (/nl)- Note that while point sources 
can introduce a significant bias on the estimate of f^™, 



they do not change significantly the error bars, A/j^^j that 
depend linearly on square root of sky coverage given the 
overall level of noise (Table FH)) . Therefore, error bars on the 
measured f-^-^^ are not mentioned further in this section, for 
simplicity. We now discuss the results obtained in Table [O] 
starting with the 100 and 143 GHz channels, dominated by 
radio sources (ij |5.2.ip and finishing with the 217 GHz chan- 
nel, where one has to account for the additional IR source 
contribution (^ I5.2.2p . 



5.2.1 100 and 1^3 GHz: effect of radio sources 

In the first two bands, 100 and 143 GHz, faint point sources 
are composed mainly of radio sources, as can be seen in 
the ERCSC. Radio sources are not clustered: their power 
spectrum is known to be flat so they act as a positive, un- 
correlated noise. To understand the effect of such a noise, we 
study in details the configuration of a minimal mask in the 
143 GHz channel as illustrated for each functional by Fig. [S] 
As a positive uncorrelated noise, radio sources do not 
affect significantly the Area MF. Indeed, the positive na- 
ture of such a noise is subtracted out when computing the 
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Table 9. Estimates of /^^ ™ presence of point sources. To create the test maps, we used the procedure described in 
eq. 12411 with a noise at the level of the extended mission in each band. The analyses are performed for combined MFs, 
V0 + V1 + V2 + V3. 



density contrast 5 of the map and, as discussed in § 3] the 
convolution effect of a zero average distribution on the Area 
is negligible as long as its variance is small compared to the 
variance of the signal. 

The presence of point sources brings an excess of posi- 
tive clusters; it also slightly shifts the curve Vk{v) to smaller 
values of 1/ as an effect of subtracting the average from the 
temperature map when computing the density contrast, as 
can be easily noticeable on the Perimeter panel of Fig. [S] 
These effects induce an overall positive bias for the Genus 
and the Perimeter, -while A'ciuster presents an excess for posi- 
tive thresholds but a deficit (of under-dense regions) for neg- 
ative thresholds. These biases remain after renormalisation, 
i.e. -when passing from Vk to Vk, except for the perimeter, 
where division by the factor At inverts the bias. Indeed, the 
presence of point sources increases the ratio ai/ao, hence 
the measured value of Ak, k > (see Appendix IX)). 

Table [TOI shows the corresponding bias on the measure- 
ment of /nl in the null hypothesis /^l™ ~ 0. The results of 
this table can partly be inferred intuitively from the exam- 
ination of Fig. |8l small bias for the Area, positive bias for 
the Genus and A'ciusters, negative bias for the Perimeter, re- 
sulting in an overall positive bias for the combination of all 
functionals. The examination of other masking levels con- 
firms the results of this analysis: the effect of point sources 
in the 100 and 143 GHz is a positive bias on the measured 
/nl which, in addition, does not depend on the primordial 
level of non Gaussianity, 

L.=f^'^ + f^\ (25) 
and decreases, as expected, when more point sources are 



^prim^O Vb Vi V2 Vs All 

(/nl) 0-8 -8 14 2.3 10 

Table 10. Bias on the measurement of /j^^ introduced by point 
sources at 143 GHz in case of weak masking at the lOo detection 
level in the ERGSC. Each column corresponds to using a specific 
functional in the analysis, or, for the last one, the combination 
of all functionals. Here, the null hypothesis /j^"™ = is tested, 
but in practice, the effective bias does not depends on the value 
of /P^i- (eq.Iig. 



excluded by the masks, as illustrated by Table 1111 In par- 
ticular, the bias induced by point sources becomes nearly 
negligible compared to expected error bars {Af^^ > 10, Ta- 
ble [T4| when masks are set up at the 3o- level. 

5.2.2 217 GHz: effect of radio and IR sources 

In the 217 GHz band, in addition to radio sources, an IR 
background contributes to the faint point sources, which re- 
sults in a new bias on the measurement of /^^ 1 Table [5] 
shows. 

The contribution from unclustered radio sources is a 
decreasing function of frequency and masking. It should act, 
as in the 100 and 143 GHz channels, as a positive bias on 
the measurement of f^^ that does not depend on the value 
of/P™ (see eq.m- 

On the other hand, IR sources are clustered and form 
mainly a diffuse, unresolved background, which cannot be 
dealt with masks. They induce a bias on the measurement 
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Figure 8. Effect of the presence of point sources on the measurement of Minkowski Functionals in the 143 GHz channel. Here, only the 
brightest point sources, with an ERCSC signal to noise ratio larger than fi = lOcr, are masked out. The configuration is the same as in 
Fig. |6l normalised functionals vy^ are plotted as functions of the threshold v, after subtracting the Gaussian limit prediction, v'^ . The 
black thick curve corresponds to f^-^ = but with point sources, while the two other one correspond to the expected curves in presence 
of primordial f^^^ as shown inside each panel. 



of /pjL which appears to depend on the value of /j^™ as 
can be inferred from Table [O] To analyse this bias more in 
depth, we concentrate on a configuration where the contri- 
bution of radio sources is masked out as much as possible 
and nearly negligible, with the 3a flux cut setting for the 
masks. Figure |5] displays the functionals obtained in two 
cases, = and = 50. It is interesting to notice 

that the curve obtained for /jSJ"™ = 50 is nearly exactly the 
sum of the curve for = 50 with no point source con- 

tribution and the curve for .f^™ ~ with point sources. 
Unfortunately, this linear property does not translate in a 
simple way in terms of bias on the measured /j^^^ , when per- 
forming the analysis, as illustrated by Tables [121 and [T3l 
What we find, instead, is the bias due to IR sources to be 



roughly proportional to /^™, when one considers the com- 
bination of all functionals to perform the measurements. 

Our final approximation for the total bias in the 
217 GHz channel is therefore 

U @217GHz = + fl^'''^'''" + fl^'''' (26) 
with jbjj's.i'^idio {jepgji(jijig only of the masking level, as in 
ij |5.2.1l Here, this bias grows from negligible for the 3a masks 
to ~ 5 for the lOa masks (right column part of 

Tabled] with ./nl™ = 0)- The other term, f^^"'^^ does not 
depend on masking but is approximately proportional to 
/■prim £ moderate values of 

/■prim 

/™^^. l/.T"i< 50. (27) 



Note that with Planck extended mission signal to noise, one 



Non Gaussianity and Minkowski Functionals 15 



0.0000 




•/S = 3(7, tNL=50 

■ /S = 3(T, W=0 

■ fM,=50 



-0.0010 L 



-2 



0.004 

! = 3<7, tNi.= 50 
0.003^ -- fi = 3a, f„,,= 




-0.003 



0.008 



3 0.004 



0.000 



-0.002 



, 1 1 1 , 

. /S = 3ff, f„,= 50 




1 , 1 1 1 , 


--/S = 3(7, f,,= 






fNL=46.5 A 






fNL=50 / \ 






/ 






/ ^ /^"^ 












/ V V 




H 




\\ 






\ \\ 




il 1 


\ \ 


'/ \ \ 


' / ' 


\ 




-^•^ / / 


\ 


Ml \\ 


/ / 


\ 


M V '** , 


V / /■ 


V 


TV y 








, , '"T 1 , 







0.004 - 



0.002 



> 0.000 



-0.006 



-2 




Figure 9. Effect of the presence of point sources on tlie measurement of Minkowski F\inctionals, similarly as in Fig. [S] but for the 
217 GHz channel. The black thick solid and dashed curves correspond respectively to primordial ^^-^ = 50 and j^-^ = 0, plus point 
sources, while the two other ones correspond to the case of primordial f^-^ only, as shown inside each panel, with f^-^ = 50 for the blue 
curves, and with the value of f^^^ found when estimating this quantity in the foreground-biased maps (Table [T3) l . in red. Note that on 
upper left panel, there is no red curve, and the black curve superposes exactly to the blue one. To emphasise the effect of the clustered 
IR background, a 3cr level masking was performed to subtract as much as possible the contribution from radio-sources. 



expects to be able to characterise more accurately the IR 
background. It might then be possible to account for it in a 
better way, by e.g. including it in the model itself instead of 
ignoring it due to lack of precise knowledge. 



6 FOREGROUNDS II: GALACTIC 

RESIDUALS AND GALACTIC MASK 

Galactic signals are a major issue for cosmological studies 
of the CMB and have to be accurately assessed. It is usu- 
ally treated in two ways: masking alone or in conjunction 
with component separation. Here, we test the effects of the 
Galactic residuals left from these methods on the estimation 
of with Minkowski Functionals by using simulations of 



the Foreground emission, masks constructed from these sim- 
ulations for different sky coverages, and a naiVe model of 
component separation quality. 

The philosophy adopted here is similar to that in the 
previous section: we do not include in the model templates of 
Galactic foregrounds. Instead, we analyse the biases on f^j^ 
introduced by neglecting the presence of galactic emission. 
Such biases are expected to be small if one restricts to re- 
gions far from the galactic plane or/and if proper component 
separation has been performed prior to the measurements. 

This section is divided into four parts. The first one, 
§ 16.11 details our simulation settings. The second one,§ 16.21 
looks at the statistical uncertainty expected on the mea- 
sured /nl from Minkowski functionals as a function of sky 
coverage. The third one, § 16.31 examines, as functions of 
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Flux cut /biiis OlOOGHz @143 GHz 
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5fT 8 4 

3(T 3 1.5 

Table 11. Bias on the measurement of f-^-^ introduced by point 
sources at 100 and 143 GHz as a function of masking level, ex- 
pressed here in terms of ERCSC signal to noise threshold, as- 
suming that the combination of all functionals is used to mea- 
sure /j^L- This bias, modelled by eq. I I25I I. does not depend on 
the actual level of primordial non Gaussianity, /JJ"™' ^'^'^ 
therefore be easily corrected for. Here, as detailed in § 15.1.31 the 
bias is obtained in the following configuration for the extended 
mission: Afgjtjc = 1024, £max = 2000, Gaussian smoothing with 
= 5', rtbins = 26 and i-'max = 3.5. With a different set up, 

c W rl ivl 

a new estimate of the bias would be needed, but this is an easy 
exercise. 
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Table 12. Point sources bias at 217 GHz for /^j^ = 0: for a flux 
cut corresponding to the 3(T detection level in the ERCSC, i.e. 
an important masking of radio sources, we see the effect of the 
clustered IR background, in the case of a null primordial non 
Gaussianity, /P™ = 0. 



sky coverage, the biases brought by Galactic foregrounds, if 
these latter would not be removed at all, and the expected 
improvements on such biases brought by component sepa- 
ration. 

6.1 Method 

To perform our simulations, -we need a -way to generate re- 
alistic maps of the Galactic emission, w hich is reviewed in 
§ 16.1.11 To do so, we use again the PSM l|Delabrouille et al] 
|2012| ). 16. 1.2 1 explains how masks are generated, in order to 
exclude regions where Galactic emission is too strong, 16. 1.31 
details the very simple method we decided to use to assess 
the affect of component separation quality. Other technical 
details of our simulations are provided in § 16.1.41 

6.1.1 Simulations of the foreground emission: the 
different components 

We use again the PSM to construct the Galactic emission. 
The code uses template maps interpolated at the desired 
frequencies. The detailed modelling of each component is 
described in iDelabrouille et al.l l|2012l '). Here are the list of 
the components we simulated in our map of the Foreground 
emission. 





Vb 


Vi 


V2 


V3 


All 


(/nl> 
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69 
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44.5 


69 



Table 13. Point sources bias at 217 GHz for /^j^ = 50: same as in 
Table [T^ but for a significant level of primordial non Gaussianity. 



The diffuse Galactic emission arises from of the inter- 
stellar medium (ISM) in the Milky Way. The ISM is com- 
posed of different phases, from cold molecular clouds to hot 
ionised regions, of magnetic fields and cosmic rays. The in- 
tensity of the corresponding emission depends on Galactic 
latitude. It is stronger near the center of the Galaxy and 
decreases at lower/higher latitudes following approximately 
a cosecant law 1/ sin(|6|). More precisely, we usually classify 
the different ISM components according to their physical 
emission processes: 

(i) synchrotron radiation is emitted by relativistic cosmic 
rays spiralling in the Galactic magnetic field. Its intensity 
depends on the cosmic ray density and on the magnetic field 
strength perpendicular to the line of sight; 

(ii) free-free emission originates from the ionised medium 
in the ISM, as a result of the interaction of free electrons 
with positively charged nuclei. It comes principally from star 
forming regions in the Galactic plane; 

(iii) there is also a thermal emission coming from dust 
grains heated by stars, which is the dominant contribution 
above 70 GHz; 

(iv) another type of emission have been found at mi- 
crowave frequencies which is probably due to small spinning 
dust particles ("Anomalous emission"); 

(v) At these frequencies, there are also molecular lines 
emerging from molecular clouds, particularly those of ^^CO 
at 100 and 217 GHz. 

For the last contribution, i.e. CO lines, templates and 
models of the emission are mostly unknown at present time 
and we choose not to model it and to simulate only the 
143 GHz part of the Galactic foreground. However, concern- 
ing noise level, our analyses correspond to a combination 
of the 100, 143 and 217 GHz channels in the framework of 
the Planck extended mission. The specific CO contribution, 
even if sub-dominant compared to the other physical pro- 
cesses, will be studied when Planck templates of CO will be 
available. 

The resulting foreground map is represented on Fig. IIOI 
with two colour scales enhancing different aspects of this 
emission. 

6.1.2 Galactic Mask 

To create Galactic masks, we consider the fraction /sky of 
the sky we aim to keep, which sets up an intensity thresh- 
old for the Galactic foreground above which the correspond- 
ing region of the sky is masked out. Once the mask func- 
tion M is set-up, which is equal to one for valid pixels and 
zero for excluded pixels, convolution of this function is per- 
formed with a Gaussian kernel of size 0f,„„.. = 5°, to ob- 

FWHM ' 

tain a smoothed version Msnioothod- New masks with smooth 
boundaries are extracted from this map, by selecting pixels 
with Afamoothcd 5S A/thrcsh and cxcludiug the others, where 
the value of Mthrcsh is tuned to get the correct value of /sky 
An example of masks constructed that way is given on Fig- 
ure [TT] 

6.1.3 Component separation quality 

In order to asses simply the impact of the residuals of 
component separation, i.e. of its quality, we model compo- 
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Figure 10. Map of Galactic emission at 143 GHz, units are in 
Kelvin. In the top panel, the colour scale is linear. In the bot- 
tom one, the colour scale is histogramequalised to increase the 
dynamic range and make visible both the regions of low and high 
emission intensity. We have included synchrotron, free-free, ther- 
mal emission and emission from spinning dust particles. 
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Figure 12. Simulations of a microwave observation at 143 GHz, 
including CMB and Galactic emission in two cases: one with no 
component separation a = 0.77 and one with a "good" compo- 
nent separation a = 0.01. Units are in Kelvin. 



itudes will be the same. The value a = 0.77 corresponds to 
the initial level of foreground emission, so it is equivalent 
to no component separation (top panel of Fig. I12[) . Real- 
istic values of a, when examining maps obtained from ac- 
tual component separation methods, rather appear to range 
typically between a = 0.01 (bottom panel of Fig. I12[l and 
Q = 0.05. 

Obviously our modelling of Galactic residuals after com- 
ponent separation is very rough, but it should be sufficient to 
assess what should be the value of /sky for making the biases 
on /pjL induced by these residuals negligible. A better mod- 
elling of the residuals would require detailed examination 
of the results obtained from actual component separation. 
Furthermore, a new analysis would be required each time a 
new component separation method is considered: this is far 
beyond the scope of this paper. 



Figure 11. Example of a galactic mask with /g^y = 0.70, drawn 
from the PSM. It relates directly to the Galactic emission shown 
in Fig. [TO] 



nent separation results by simply adding to the CMB the 
map of Galactic foregrounds multiplied by a scaling factor 

C"CMB(QHL , J i. 1 

a, where acMBiaHL and crpGaHL are respectively 

CpGiaHL 

the standard deviations of the CMB and foreground maps at 
high latitudes, i.e. measured in pixels outside the mask with 
/sky = 0.80. So for a "quality factor a — 1, the normalised 
contributions of Galactic foreground and CMB at high lat- 



6.1.4 Simulations pipeline 

The simulation strategy used here is the same as in § 15.1.31 
the Galactic foregrounds are added only to the test maps, 
or, in other words, the "data" maps. The rriNG = 200 maps 
used to compute the model prediction and the m — 10 000 
used to estimate the covariance matrix in the function do 
not have the foregrounds, but the treatment is the same oth- 
erwise. We test 4 values of /p™ = {-10, 0, 10, 50}, for each 
of which mtest = 200 test maps are generated, to which we 
add Galactic foregrounds in the / — 143 GHz channel, with 
a level a as described in § 16.1.31 Then Gaussian beaming is 
performed with the beam of the instrument in this channel. 
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Table 14. Sensitivity of the estimator versus sky coverage: 
(Af )^/^ is the error bar estimated from the combination of 
all functionals, Vq + Vi + V2 + Vz- It is approximately a linear 
function of \//sky • The estimates of the error on /nl a.re per- 
formed for the combination of the 3 channels in the extended 
mission configuration. 
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Table 16. Galactic bias at high latitudes, fGal-high ^ ^ f^j^,, 
tion of sky coverage, /sky, as explained in the main text. The 
results displayed on this table assume Gaussian smoothing with 
6^ = 5' in the simulations used to perform the calculations 
(§[6X1. 



the Galactic foregrounds, while § 16.3.31 examines the biases 
as functions of component separation quality factor, q. 



^FWHM ~ T-"^ 1 that corresponds as well as to the beam for 
the 3 combined channels (that we assume to be used for this 
analysis) and the noise at the level of the extended mission 
(for the 3 combined channels) is added (see Appendix ID)|. 
Then filtering is performed using either a Gaussian window 
of size or Wiener filters discussed in ij |3.3l Finally, the 

Galactic mask calculated as in S 16.1.21 is added, for a given 
/sky. Note that at variance with the point-source analysis, 
the Galactic mask is not inpainted, due to its rather large 
size: the pixels inside the mask are just ignored by the MPs 
code. The complete pipeline is summarised as follows: 



map = CMB(/P™)*beam[0^„^„(/)1 

+ a X ^2M22IiIi X Gal. FG(/) * beam[e^ 

CFGiBHL 

-f noise(/) 

+ smoothing [Sp^jj^] or Wiener filtering 



4- galactic mask (/. 



sky J 



,(/)] 



(28) 



6.2 Sensitivity versus sky coverage 

First, we test the sensitivity of our estimator to sky cover- 
age, using maps of the CMB without Galactic foreground i.e. 
a = and just looking at the resulting error bars, A/^^ . Ta- 

ble ll4l shows (A/^^)^/'^ as a function of /sky, obtained from 
the combination of all four functionals for various Gaussian 
smoothing scales ^^^hm different Wiener filterings. As 
expected, (A/j^^^)^^^ is approximately a linear function of 
^ /sky . Note that these results do not change significantly in 
the presence of Galactic foregrounds as long as /sky < 0.80. 



6.3 Effect of Galactic foreground 

In the following, we study the actual effects of the Galac- 
tic foregrounds, which present a rather complex behaviour 
as a function of sky coverage, /sky Our analyses restrict to 
Gaussian smoothing and put aside Wiener filtering. From 
the qualitative point of view, the latter indeed leads to re- 
sults similar to what is obtained with the smallest Gaussian 
smoothing scale. This section is divided into two parts. Sec- 
tion 16.3.11 assumes no component separation and analyses 
in details the biases brought on the measurement of /^^ by 



6.3.1 Two behaviours, two components 

We now examine what kind of biases the Galactic fore- 
grounds induce on f-^-^. Table [15] shows our forecast for 
(./nl) ^ * function of /sky and f^^l™^, at the smallest Gaus- 
sian smoothing scale, S^,„„.. = 5'. Remind that Galactic 
foregrounds are completely present, without any attempt to 
remove them with component separation. 

Numerous interesting results can be extracted from Ta- 
ble [HI 

(a) First, close to the galactic plane, so for /sky close to 
unity, a very important signal, that we denote by /j^^'"^'^"'^' 
overrides the primordial one. Indeed, when /sky = 0.90, 
(/nl) does not depend much on the value of /^"™ and can 
be approximated as follows: 



(./nl) = 



Gal — pla 
NL 



/■prim 

/-Gal — plane, const , J nl 

- /nl + 



with y*^^^~p^^^'-'"'-'-'^^^^ 



/sky = 0.9 

(29) 



-22 for the Vi + V2 + Vz combina- 
tion and /Gai-pianc,const ^ _26 for the VQ + V1 + V2 + Vz 
combination. 

(b) At higher latitudes, a second type of signal less power- 
ful than /^/''"P'^"° appears. This new contribution, denoted 
by /nx'"*^'^ (/sky), brings a positive bias on the measure- 
ment of /jjL s-iid does not hide the primordial signal. For 
/sky = 0.8, both /°j^'''"P'''°'= -which brings a negative bias- 
and /^jj''"'"^'' -which brings a positive bias- contribute: 



(./nl) = ./, 



Gal — plane 



Gal — high^^ Q-^ , ^prim r ri i 

NL ' ./sky - U.. 



^'^•^(0.8) + /^r 

(30) 

with /^.^'-'^'^'^(/sky = 0.8) ~ 20 and /^^''^'^"^ given by 
eq. (|29p . Those two signals compensate each other and it is 
remarkable to see that if only three functionals are used, ex- 
cluding the Area, the bias on f-^-^^ is almost negUgible (apart 
from the /nl™/10 contribution)! 

(c) Finally, for smaller /sky, the signal from the Galactic 
plane is totally hidden and only /^^'"'^'^''(/sky) contributes 
as a positive bias: 

(/NL)=/NL'-''^'(/eky) + /;^r. /sky < 0.7. (31) 

The quantity /^^'"'''^''(/sky) is shown in Table[l6]for various 
sky coverages. As expected, /^L'~'^'^'^(/sky) is a decreasing 
function of /sky 
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Table 15. Effect of Galactic foregrounds on the measurement of /nl as a function of 
sky coverage /^jjy and level of primordial non Gaussianity /^"™. The results displayed 
on this table assume Gaussian smoothing with 8^ = 5' in the simulations used to 

" FWHM 

perform the calculations (^ 16.1.41 1. 



6. 3. 2 Smoothing 

To understand more deeply the effects of the two types of 
Galactic foregrounds we found above, we examine what hap- 
pens when the smoothing scale is varied. Table flTl gives {/^^ ) 
as a function of ,„„.. , for various sky coverages. Obvi- 

FWHM ' ■.' ^ 

ously, one can infer from this table that the positive bias 
increases with smoothing scale, suggesting that the Galac- 
tic signal at high latitude dominates at large scales, while 
the Galactic plane signal is present with its negative bias 
only at small scales. This result, in addition to the analy- 
ses performed in previous paragraphs, show that Minkowski 
Functionals remain helpful in understanding and isolating 
the different biases induced by Galactic foregrounds, simi- 
larly as for the bi-spectrum. This demonstrates again the 
discriminative power of MFs. 



6.3.3 Component separation 

The analyses in previous paragraph were performed in the 
most pessimistic case, when no component separation is per- 
formed, i.e. a = 0.77 in eq. (|28p . Now we consider the 
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Table 17. Galactic foreground bias as a function of smoothing 
scale for various sky coverages, /^kyi in the null hypoth- 

esis /^"™ = 0. The measured quantity {/^l) is given for the 
combination Vq 4 Vi 4 V2 4 Vs , 



case when the Galactic foregrounds have been largely re- 
moved with a component separation method and that only 
a fraction remains, a < 0.77. Table [18] shows the bias ex- 
pected on the measured fj^^^ due to Galactic foreground 
residues as a function of a and for various sky coverages. 
Here, we focus on small scales, with Gaussian smoothing 
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Table 18. Galactic foregrounds bias as a function of component 
separation quality factor a and sky coverage. The forecasted value 
(/nl) given in the null hypothesis /JJ"™ = 0. Note that this 
bias would also stand approximately for a value to add to Z^"™ 
in the case /P"™ 7^ 0. 



at S^,.,„., = 5' and on the measurement of f.,, with the 

FWHM 

combination Vo + Vi + V2 + V3 . 

Table [TH] shows that with a good but realistic compo- 
nent separation (0.01 < a < 0.05), the bias due to the Galac- 
tic foregrounds becomes negligible compared to error bars 
(Table [HI) even for rather large sky coverage, /sky = 0.8. 



7 CONCLUSION 

In this article we have studied in detail the ability of 
Minkowski Functional^ of excursions of the temperature 
fields to estimate primordial non Gaussianity, f-^-^, in a 
Planck like experiment. To do that we used a standard 
Monte-Carlo approach to define a statistics assuming 
the weakly non Gaussian regime, where the errors are dom- 
inated by the Gaussian part of the signal. We first assessed 
the numerical limits in the approach. Then we studied 
the effects of inhomogeneous noise, point sources and galac- 
tic foregrounds on the determination of /^^ ■ The main re- 
sults of our investigations, performed for the 100, 143 and 
217 GHz channels, are the following: 

(i) It is best to measure normalised functionals, Vk = 
Vfc/Afc, to reduce the effects of inhomogeneous noise. 

(ii) The functionals are all sensitive to non Gaussianity, 
but the following hierarchy can be set: Perimeter > Genus 
> A'^ciustors ^ Area. It is worth combining several Minkowski 
functionals to obtain better constraints on f-^-^ , although the 
Area does not improve much the results. 

(iii) To extract most of the information of interest while 
keeping the approach valid, it is convenient to perform 
the analyses in the "3.5 sigma excursion range with a num- 
ber of bins ribins ~ 26 for each of the functionals. 

(iv) To have proper convergence of the x^ , it is required in 
practice to perform about m = 10 000 Gaussian simulations 
to estimate properly the covariance matrix. 

(v) Combining Wiener filtering for the field and its deriva- 
tive bring the best constraints on f^^^ when using MPs. 
In particular, Wiener filtering does better than Gaussian 
smoothing. Note that most of the information on non Gaus- 
sianity is contained at the smallest possible scales, of the 
order of the beam size. 



^ and the number of clusters, that we call a Minkowski Functional 
here to simplify the presentation. 



(vi) Point sources foregrounds introduce a bias on the 
measurement of /^^l that can be estimated and thus cor- 
rected for accurately. Note that with appropriate masking 
of the brightest sources followed by a simple in-painting pro- 
cedure, this bias becomes negligible except at the 217 GHz 
frequency. 

(vii) Galactic foregrounds introduce a complex bias on 
the measurement of f^^ that depends on the fraction of 
sky covered, /sky, or in other words, depends on how much 
the most luminous part of the galaxy has been masked out. 
This bias can be again corrected for and is coincidentally 
negligible when the combination of Perimeter, Genus and 
number of clusters is used and /sky ~ 0.8. With appropriate 
component separation the bias due to Galactic foregrounds 
should become negligible compared to the error bars on /^^ > 
at least for /sky < 0.8. 

(viii) With all the effects described above under control, 
we expect to be able to measure /^^ using Minkowski Func- 



tionals with an error of the order of A/„ 



10. 



These are excellent news overall, since it suggests that 
Minkowski Functionals should be capable of putting inter- 
esting constraints on /^^ > even in view of the performance of 
fully optimal estimators for that purpose. It is reasonable to 
assume that they should also provide non-trivial constraints 
on other types of primordial non-Gaussianity even if there 
are no prediction yet at the time of the analysis which would 
allow building dedicated and somewhat more stringent indi- 
cators. Finally, although rather detailed, our analyses could 
be improved by lifting the following limitations before prac- 
tical applications are considered: 

(i) Our model for testing component separation quality 
is rather naive and two major issues remain to be addressed 
before drawing definite conclusions: (a) component separa- 
tion methods do not remove galactic components the same 
way in each part of the sky nor each component: the descrip- 
tion of Galactic residuals in terms of a contribution simply 
proportional to them has to be improved; (b) here, galac- 
tic components are "added" everywhere: component separa- 
tion can subtract CMB signal too, an effect that we did not 
consider here. In fact the analysis of component separation 
met hod quality needs a specific study for each method at 
use iLeach et al.l |200^). 

(ii) In this study we did not co nsider the use of foreground 
templates in the reference maps (|Komatsu et al.|[2002l [2OI1I ) 
to assess for the residuals. The use of templates, if they are 
reliable, is expected to correct for the biases introduced by 
the foregrounds. We can thus consider the results of our 
analyses as pessimistic in that respect. 

(iii) We did not characterise biases induced by sec- 
ondary anisotropics: they could be important, even dom- 
inant, as a bi-spectrum is created from the covari- 
ance between weak lensing and Sunyaev-Zeldovich effect 
or Integ rated Sachs- Wolfe effect ( I SW). Indeed, previous 
studies JColdberg fc Spergell 1 19991 : ISerra fc Cooravl I2OO8I : 
iHanson et al.ll2009l ) have warned about these spurious sig- 
nals and a future work about their effects on Minkowski 
functionals is planned. 

(iv) Our analyses neglected correlations in the noise. In a 
future work we shall refine them when reliable simulations 
of correlated noise are available. 
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In this paper, we have not used explicitly the analytic for- 
mulations of Minkowski functionals for our n on Gaussian- 
ity s tudies, as what was done for example in iHikage et alJ 
1 2006. 2008). Ac tually, using the "Skewness parameters" 
l Matsubarall2003l ') to study biases induced by secondaries, 
galactic foreground and point sources would allow us to use 
the numerous bi-spectra studies in an effective way. On the 
other hand, the advantage of our Monte-Carlo method is 
that it can be generalised to any stati stics, e.g. for instanc e 
the skeleton length in the excursion (|Novikov et al.ll2006l ). 
or any type of n on Gaussianity, e.g. t hat induced by cos- 
mic strings (e.g., iBouchet et al. 1 1 19881 : iRingevall |2010| . and 
references therein). 
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APPENDIX A: MINKOWSKI FUNCTIONALS: 
DEFINITIONS AND THEORY 

For a field f{x) of zero average and variance ctq defined on 
the two-dimensional sphere S^, an overdense excursion set 
writes 

E = {t G J{x) > uao}. (Al) 
The boundary of the excursion is 

dT. = {xe S^l fix) = uao}. (A2) 
Then the three Minkowski functionals on the sphere write 

(A3) 

dl, (A4) 



Area : Vo{v) = — / dQ. 

47!- Je 

Perimeter : Vi(u) = / 

47r 4 Jgj 



(A9) 



The amphtude Ak depends only on the shape of the power 
spectrum Ct. 



UJ2 



2 



k<2 



(AlO) 



(All) 



where cok = ■k''^'^ /V[k/2 + 1), which gives ujq = 1, = 2, 
CJ2 = TT and GO and o\ are respectively the rms of the field 
and its first derivatives. 

In the weakly non Gaussian regim e, one can write, at 
leading order (see, e.g. lMatsubarall2O10h 



(0) , (1) 



(A12) 



where w^"-' is given by eqs. (|A7|I and (|A8|I . while the first 



order non Gaussian correction writes 



kSi 



k{k - l)Su 



(A13) 



with 

and where 



' dv = e ' \ — erfc 



x/2 



S = 



Si 



Sii 



if) 



sdv/pv^/) 



(A14) 

(A15) 
(A16) 
(A17) 



Genus : ¥2(1') — 



_L_L 

47r 2n 



ndL 



as 



(A5) 



where dQ and dl are respectively elements of solid angles 
(surface) and of angle (distance), k is the geodesic curvature. 
Note that the Genus can be also expressed as the number of 
component^ in the excursion minus the number of holes in 
the excursion. 

The fourth functional we use in this paper, ¥3(1'), is 
defined, for v > 0, as the number of components in the 
excursion. Symmetrically, for < 0, it is the number of 
underdense components (or the number of components in 
the excursion {x G f{x) < uao}). 

In the Gaussian limit, th e fu nctionals can be expressed 
the f ollowing way (see, e.g. iMatsubara .2010. : .Vanmarcke 
1 19831 ): 



with 

V3{u) = 



Vk()y) = AkVk{)y), 
exp(-//V2)-fffc-i(!^), k^2 



erfc {iy/V2) 



(A6) 

(A7) 
(A8) 



A component is a connected subset of the excursion. 



are skewness parameters (|Matsubarall200 . Each of them is 
a weighted integral of the bi-spectrum and is thus directly 
proportional to /nl. 



APPENDIX B: ALGORITHM USED TO 
COMPUTE MINKOWSKI FUNCTIONALS 

The code we developed is available by simple e-mail request. 
The numerical technique we use for measuring Minkowski 
functionals on HEALPix maps, consists, for each threshold 
value of the temperature, of two steps. 

In the first step of our algorithm, an on-grid cluster anal- 
ysis is performed in order to define the connected ensembles 
of pixels that have temperature values above the threshold. 
The output of this step is a map of integers (flags) that as- 
sign a negative "outside" flag for pixels below the threshold 
and a positive flag that corresponds to the cluster number 
of each connected region above the threshold. Optionally, 
the boundary pixels (defined either within high excursion 
regions or just outside them) can be marked. The implemen- 
tation follows closely the Cn d_reg2d procedure developed 
first f or the Adhesion model J Poeos vani 1 1 9891 : iKofman et al.l 
1 19901 . 1 19921 ). 

The on-grid cluster analysis then proceeds as follows. 
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At the onset, all pixels are considered as "unseen". In the 
outer loop, the code scans the map, checking the temper- 
ature values until the first pixel exceeding the threshold is 
encountered. Along the way, the checked pixels that have 
been found to be below the threshold are marked with the 
"seen, below" flag. The new pixel above the threshold is as- 
signed the flag that corresponds to the cluster order number 
(zero for the first found). Its neighbours are investigated. 
The ones below the threshold are marked as "seen, below" , 
the ones above the threshold are given the same cluster num- 
ber flag and are put onto stack for further analysis. Next the 
pixel from a stack is drawn, its neighbours are checked in 
the same way, with the ones above the threshold further 
added to the stack. This inner loop proceeds until the stack 
is exhausted which signifies that all connected pixels belong- 
ing to the first cluster are found. The control is reverted to 
the outer scanner that proceeds checking the pixels, skipping 
over those already analysed, until the next previously unseen 
pixel above the threshold is found. This pixel acquires the 
next cluster order number and then the inner, stack-based, 
loop finds all the pixels connected to it, and so on. 

This code is fast and linear in the number of pixels since 
the neighbours of each pixel are analysed only once. Option- 
ally, the boundary pixels can be marked differently from 
pixels inner to the regions. Different decisions about what 
constitute connected pixels can be easily implemented. In 
the current implementation, we consider all pixels that have 
at least a common vertex belonging to the same cluster (as 
we shall see below this has to be taken into account when 
considering the Euler characteristic of the total excursion 
set). The immediate by-products of the cluster analysis are 
the volume (area) of all connected regions above the thresh- 
old, which constitutes the first Minkowski Functional, and 
their number iVciustcrH 

The other two Minkowski functionals of the 2D excur- 
sion sets, the Euler characteristics and the perimeter of the 
set are computed in the second step of our algorithm, using 
the just obtained clustering information. 

The Euler characteristic of each individual cluster, and, 
after summation over all the clusters, of the whole excur- 
sion set above the threshold, can be obtained in one pass 
by analysing the pixel grid vertices on the cluster bound- 
aries. Gauss-Bonnet theorem links Euler characteristics of 
the region to the integration of the curvature of its bound- 
ary. However, since topological properties are invariant un- 
der any continuous transformation of the boundary, the need 
to explicitly evaluate the boundary curvature is eliminated 
by considering the curvature to be accumulated just in the 
outside vertices of the boundary pixels. Thus one only needs 
to assign the appropriate curvature weights to the grid ver- 
tices and sum over them. Si milar procedure o n the Cartesian 
grid has been described in lGav et al.l l|2012l ). A 3D version 
of the code is also readily available. 

Clearly, the vertex contribution is determined solely by 
the temperature distribution of the pixels that form this 
vertex - which are below and which are above the threshold. 
Necessary weights can be boot-strapped by considering the 
elementary situations. Most of the grid vertices in HEALPix 



° To compute A^ciustcr for negative thresholds, one just multiplies 
the map by -1, and repeats the procedure. 



pixelisation are regular, being formed by four adjacent pix- 
els. 

(a) Consider a single pixel cluster above the threshold. 
Its Euler characteristic is x = 1- It has 4 boundary vertices, 
all four equivalent, formed by one adjacent pixel above the 
threshold and three below. They should contribute equally, 
thus we assign the weight 1/4 to any vertex of this type. 

(b) Consider next a two pixels cluster, formed by pix- 
els having a common side. The boundary has six vertices, 
four corner ones of the type (a) and two new ones, which 
are formed by two side-by-side adjacent pixels above the 
threshold and two below it. The cluster has x = !> thus new 
vertices must contribute a weight equal to 0. 

(c) Consider a three pixels cluster that forms a corner. It 
has 5 vertices of type (a), which weights add up to 5/4, two 
vertices of type (b) and one new vertex that is formed by 
three pixels above the threshold and one below it. For the 
total X = 1 this vertex must contribute the weight —1/4. 

(d) The last possibility arises when we consider the 
boundary vertex formed by two pixels above the thresh- 
old that are touching just at the vertex. Its weight depends 
whether we consider the clusters to be linked through the 
vertex or disjoint. In the former case, which corresponds to 
our choice in cluster analysis, the weight is 1 — 6/4 = —1/2, 
while in the latter case it is 2 — 6/4 = -1-1/2. However, for 
statistical analysis of the total Euler characteristic of the ex- 
cursion sets, both choices, while exact for the corresponding 
clustering decision, will lead to biased results. Indeed, one 
may argue that in a situation when the discretised field has 
high and low pixels mixing at a vertex, it is equally proba- 
ble that the regions below the threshold or above it connect 
through this vertex. Assigning the weight zero in this case 
will reflect such a "statistical" consideration and this is the 
choice that we make in our statistical analyses. 

(e) The vertex that is formed by four pixels above the 
threshold is not a boundary vertex and its weight is zero. 

The HEALPix grid on a sphere has eight special vertices 
that are formed only by three pixels. Their weights need 
to be defined separately. Consider elementary clusters that 
have one of such vertices. 

(f) A one pixel cluster with a special vertex has three 
regular vertices of type (a). Thus, a special vertex with one 
pixel above the threshold and two below has the weight 1 /4. 

(g) A two pixels cluster with a special vertex on a side 
shows that a special vertex formed by two pixels above the 
threshold and one below has the weight 0. 

(h) A three pixels cluster with a special vertex inside has 
only 3 exterior regular vertices of type (a). Thus, a special 
vertex formed by all three pixels above the threshold, albeit 
being an inner one, has a weight 1/4. These eight special 
vertices are providing the Euler characteristic of the entire 
spherical manifold, x = 2, when all pixel values lie above 
the threshold. 

Once the weights are defined, the computation of the 
Euler characteristic of individual clusters and the entire ex- 
cursion set is a simple one pass loop over vertices to deter- 
mine their type from the temperature value at the adjacent 
pixels and add the weights based on the index of pixels that 
are above the threshold. Our method does not involve any 
differentiation nor integration of the field. Moreover, if only 
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the total X needed, the cluster analysis step can be omit- 
ted. The code works in the presence of masks of arbitrary 
complexity, treating the field as defined on a manifold that 
itself has a non-trivial Euler characteristic. 

The measurement of the perimeter of the excursion set 
is somewhat more complex. Our technique is again based 
on a scan over the boundary grid vertices, for each of which 
we linearly interpolate the field based on the values in the 
four (or three) adjacent pixels. Interpolation is least-squares, 
linear, performed in the plane tangent to the sphere at the 
vertex. Interpolation requires first derivatives of the field 
which are precomputed by standard HEALPix routines and 
stored. After finding the best-fit plane for the local field 
behaviour, we compute the intersection of this plane with 
the polygon formed by the pixel centres around the vertex. 
The length of the intersection is then added to the perimeter 
of the cluster the vertex belongs to. We found that linear 
interpolation is sufficient given the level of accuracy offered 
by the pixelisation and there is no advantage to go for a 
higher order interpolation in comparison with increasing the 
pixelisation level. 

The code we developed along these lines is available by 
simple e-mail request. 



APPENDIX C: INVERSION OF THE 
COVARIANCE MATRIX: A CONVERGENCE 
STUDY OF THE 

The inversion of matrix C in eq. ((Tjl is performed using Sin- 
gular Value Decomposition (SVD): C — UDU^ or 



Cij — ^ ^ UikdkU. 



jk , 



(CI) 



hence 



=1 j=i fe=i 



(C2) 

where U is an orthogonal matrix and D a diagonal ma- 
trix such that its diagonal terms dt, the singular values, are 
ranked in descending order. Here, these latter are expected 
to be positive, since C is definite positive by construction. 
When there is a large contrast between the singular val- 
ues, r = di/dn ^ I, the calculation of the inverse of C 
becomes an ill-conditioned problem if C is not computed 
with sufficient accuracy. In our calculations, the accuracy is 
mainly controlled by the number m of simulations used to 
estimate C, but it also depends on other parameters such 
as the number of bins, the value of I'max, the smoothing 
scale £, the presence of masks, etc. In particular, the ratio 
r increases when the space between successive bins in the 
excursion is reduced, because they become more correlated. 
A sufficiently accurate calculation of C can be costly. For in- 
stance, it takes about 3 mn at present time to mesure MFs on 
a HEALPix map with A^sidc = 2048, i^max = 3.5, T^bins = 26 
and a Gaussian smoothing of Sp^m, = fO', an operation 
that can become prohibitive if repeated many times, this is 
why it is worth investigating in details what should a reason- 
able value of m to have a few percent error on the estimate 
of the x^- 

In this appendix, instead of studying the accuracy in 



the determination of the inverse of C, we analyse directly 
the convergence of the as a function of m for /nl ~ 0. 
We focus on the perimeter functional, y = Oi, but simi- 
lar results would be obtained for other functionals. To per- 
form the Gaussian simulations for computing the covari- 
ance matrix, we shall use the following typical configura- 
tion: the HEALPix resolution is A'^sidc = 2048 and the 
power-spectrum is calculated as detailed in Appendix [D] 
The maps are convolved with a Gaussian beam of size 
^FWHM ~ ) then supplemented with a white noise of 
0.33/iK.deg, and finally smoothed with a Gaussian smooth- 
ing kernel with Sp„m, ~ fO'. Our concept of infinity corre- 
sponds to m = 20 000 simulations for estimating the covari- 
ance matrix in the "ideal" case. 

Here we want to be able to estimate m for various 
choices of the parameters which are the most infiuent in 
the convergence of the covariance matrix, namely the num- 
ber of bins, n and the excursion range, J^max. We already 
computed a covariance matrix C20*: for a rather large value 
of m = 20 000. However, the optimal value of m remains 
unknown and might be larger than 20 000 in some cases. 
Now, to make a convergence study, we should realise a large 
number of trials, each of them with a subset of m reali- 
sations, in order to estimate the error in the estimate of 
C with m simulations. To avoid that procedure to be pro- 
hibitive, instead of generating random maps and measuring 
the genus on each of them, we use the property that the ran- 
dom vectors x = y'^ — y'^ are nearly Gaussianly distributed 
to simulate directly and rapidly a Gaussian distribution of 
m vectors x of zero average and covariance C20fe. To do that, 
we use the following standard procedure. We draw 



An = N{Q, I) of size n, 



then the vector 



X = U{y/D.x^) 



n n n 



i£^) = UDU"^ = Caofe. 



(C3) 
(C4) 

(C5) 



The point is to see whether using m realisations, x 
I, • ■ • , m, of vector x to estimate the covariance matrix is 
sufficient to compute the x^ accurately enough. For one data 
realisation with /^^ = Oi 



- _ -(0) , -G 
y = x'-'+y , 



(C6) 

where i'"^ is generated randomly the same way as just 
above, the estimator of the writes 

XMcim) = [y- y{Ui^)f CMc("^) [y - vUni.)] . (C7) 
with the covariance matrix estimated from m realisations. 



C:Mc(m).l^i«[f«]- 



(C8) 



while the "true" x^ should write 

xLc = [y - y{f.,)f C-'[y~ y(/,J] , (C9) 

where we approximate here by the inverse of Caofc. The 
relative difference between the exact x^ and the one driven 



from m simulations writes 



Xmc {m) - x?r 

Xtruc 



(CIO) 
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where x = {x'^'^' , x*"^' , ■ • ■ , a;'™'). To assess the convergence, 
we need the second moment of over many reahsations M 
of X [this corresponds to M x (m + 1) x n random numbers 
in total] to be smaU 

^e'. (CU) 

We performed this exercise for e = 2%, with M = 1000 and 
for the following set of values of m = 1000 i, i £ {1, . . . , 20}. 
The results of our analyses are summarised in Table [1] 

APPENDIX D: SIMULATING PLANCK DATA 
Dl Observations with Planck 

Planck is the third and latest generation of space obser- 
vatory designed to observe the anisotropies of the cosmic 
microwave background (CMB) over the entire sky, at a high 
sensitivity and angular resolution. It observes the sky in 
nine frequency bands: the Low Frequency Instrument (LFI; 
iMandolesi et al.ll20ld ) covers the 30, 44 and 70 GHz bands 
with amplifi ers cooled to 20 K; the High Freque ncy Instru- 
ment fHFI: iPlanck HFI Core Team et ai]|201ll ') covers the 
100, 143, 217, 353, 545 and 857 GHz bands with bolometers 
cooled to 0.1 K. 

In order to assess characteristics of Planck which are 
particularly relevant for the specific purposes of this paper, 
i.e. the noise level and the resolution of the CMB map anal- 
ysed, we consider the expected performances of the three 
HFI channels at 100, 143 and 217 GHz, and we adopt two 
durations for the mission: 15 months and 30 months which 
we shall refer to in the following as the nominal and extended 
mission. 

The specifications of the three channels that are studied 
here are gathered in table lDl] The beams are assumed to be 
isotropic ("circular) Gaussian. The beam size in table [DTI is 
given in terms of "FWHM" scale as detailed in § ID2I 



that paper use A'aide = 1024 and assume a truncation of the 
harmonic modes a.t £ — i'max = 2000, but we also examined 
other values of A^sidc and £max as shown on Table 21 

For the non Gaussian part of the CMB, we used updated 
versions of simulations o f harmonic coefficients provided by 
lElsner fc Wandeltl (|2009l ) which allow us to choose Nsi^c and 
^max up to ^ = 3500. 

Convolution, in particular with the Gaussian window, 
was performed in harmonic space using HEALPix package. 
In our conventions, the Gaussian kernel, writes, in harmonic 
space, 

Wi = eM-£{£ + 1)^] (D2) 

where 9s is the Gaussian size. Here we will use instead 
the corresponding full width at half maximum (FWHM): 



D2 Simulations 

The simulation procedure used in this article can be sum- 
marised as follows 

map = CMB(af,m, /nl) * beam[6lFWHM(/)] 

+ Galactic Foreground * beam[0p\iVHM(/)l 

+ Sources Foreground * beam[5FWHM(/)] 

+ noise(/) -I- Point Sources Mask inpainted(/3, /) 

+ smoothing(6|wHM) or Wiener filtering 

+ galactic Mask(/sky), (Dl) 

where / is the considered frequency. We mainly describe 
here how the first line of this equation is dealt with, while 
other aspects are examined in various sections of the article. 

To generate the Gaussian part of the CMB, we use har- 
monic coefficients a^,™ derived from standard Cold Dark 
Matter (CDM) cosmology, with th e best parameters ob - 
tained from WMAP7+BAO+i?o (|Komatsu et al.l bOllT ): 
Qa = 0.728, S^cft^ = 0.1123, fib/i^ = 0.0226, Ho = 70.4 
km/s/Mpc, Us = 0.967, r = 0.085 and A^(fco) = 2.42 x 
10-^ 

The maps are generated using HEALPix pixelisation 
jGorski et al.ll2005l ). Most of the calculations performed in 



Ducout, F. R. Bouchet, S. Colomhi, D. Pogosyan and S. Prunet 



Channels 


Beam size 9^^^^^ 


Noise (nominal mission) 


Noise (extended mission) 


100 GHz 


10' 


l.l^tK.deg 


0.7AtK.deg 


143 GHz 


7.2' 


0.7AtK.deg 


O.S^K.deg 


217 GHz 


5' 


1.1 /xK.deg 


0.7AtK.deg 


100+143+217 GHz 


7 2' a_ 


0.5 fiK.deg 


0.33 fiK.deg 



Table Dl. Planck characteristics, in terms of resolution and levels of noise, used in this 
study. 



The beam for this combination of channels is not known but should be in-between those of 
the 143 GHz and 217 GHz channel, which are the most CMB constraining at high and very 
high angular resolution. We have conservatively used the resolution of the 143 GHz for this 
study 



